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Abstract 

In this paper, we consider the stability analysis of large-scale distributed net¬ 
worked control systems with random communication delays between linearly 
interconnected subsystems. The stability analysis is performed in the Markov 
jump linear system framework. There have been considerable researches on 
stability analysis of Markov jump systems, however, these methods are not 
applicable to large-scale systems because large numbers of subsystems result 
in an extremely large number of the switching modes. To avoid this seal- 
ability issue, we propose a new reduced mode model for stability analysis, 
which is computationally efficient. We also consider the case in which the 
transition probabilities for the Markov jump process contain uncertainties. 
We provide a new method that estimates bounds for uncertain Markov tran¬ 
sition probability matrix to guarantee the system stability. The efficiency 
and the usefulness of the proposed methods are verified through examples. 

Keywords: Large-scale distributed networked control system, Markov 
jump linear system, switched system, random communication delay 


1. Introduction 

A networked control system (NCS) is a system that is controlled over 
a communication network. Recently, NCSs have attracted considerable re¬ 
search interests due to emerging distributed control applications. For exam¬ 
ple, the NCSs are broadly used in applications including traffic monitoring, 
networked autonomous mobile agents, chemical plants, sensor networks and 
distributed software systems in cloud computing architectures. Due to the 
communication network between subsystems, communication delays or com¬ 
munication losses may occur, resulting in performance degradation or even 
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instability. Therefore, it has led various researches to analyze the NCSs with 
communication delays hi, h, a, a, 0 , 0 , ra, 0 . In particular, [6] con¬ 
structed a switched system structure for the analysis of NCS by including 
actuators, sensors, and the plant as a single system. 

In this paper, we study distributed networked control systems (DNCS) 
with a large number of spatially distributed linear subsystems (or agents). 
For such large-scale systems, our primary goal is to analyze system stability 
when random communication delays exist between subsystems. Typically, 
such delays have been modeled as Markov jump linear system (MJLS) |6j, 
0, m, HU, H2, H3, in which switching sequence is governed by a Marko¬ 
vian process. Therefore, stability analysis in the existence of communication 
delays has been performed in the MJLS framework HU, Ha, |16 j. However, 
these results are applicable to the systems with a small number of switching 
modes a, Ha, Ha, [10], whereas the large-scale DNCSs in which we are 
particularly interested give rise to an extremely large number of switching 
modes. For such systems, previous approaches for stability analysis are com¬ 
putationally intractable. Although HU investigated the massively parallel 
asynchronous numerical algorithm by employing the switched linear system 
framework that circumvents the scalability issue with respect to the large 
number of the switching modes, it is developed for the independent and 
identically distributed (i.i.d.) switching. In addition, we are also interested 
in systems where the transition probabilities are inaccurately known as in 
[T6| . pL8] . [19] because, in practice, it is difficult to accurately estimate the 
Markov transition probability matrix that models the random communica¬ 
tion delays. 

This paper provides two key contributions to analyze the stability of the 
large-scale DNCS with random communication delays. Firstly, we guarantee 
the mean square stability of such systems by introducing a reduced mode 
model. We prove that the mean square stability for individual switched 
system implies a necessary and sufficient stability condition for the entire 
DNCS. This drastically reduces the number of modes necessary for analysis. 
Secondly, we present a new method to estimate the bound for uncertain 
Markov transition probability matrix for which stability is guaranteed. These 
results enable us to analyze large-scale systems in a computationally tractable 
manner. 

Rest of this paper is organized as follows. We introduce the problem for 
the large-scale DNCS in section 2. Section 3 presents the switched system 
framework for the stability analysis with communication delays. In Section 
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4, we propose the reduced mode model to efficiently analyze stability. Sec¬ 
tion 5 quantifies the stability region and bound for uncertain Markov transi¬ 
tion probability matrix. This is followed by the application of the proposed 
method to an example system in section 6, and we conclude the paper with 
section 7. 

Notation: The set of real number is denoted by M. The symbols || • || 
and || • ||oo stand for the Euclidean and infinity norm, respectively. Moreover, 
the symbol #(•) denotes the cardinality - the total number of elements in 
the given set. Finally, the symbols tr(-), p(-), <8), and diag(-) represent trace 
operator, spectral radius, Kronecker product, and block diagonal matrix op¬ 
erator, respectively. 


2. Problem Formulation 


2.1. Distributed networked control system with no delays 
Consider a DNCS with discrete-time dynamics, given by: 

Xi(k + 1) = ^ AijXj(k), i = 1,2,... ,N, (1) 

jeM'i 

where A; is a discrete-time index, N is the total number of agents (subsys¬ 
tems), Xi G K n is a state for the i th agent, M% is a set of neighbors for Xi 
including the agent Xi itself, and A,j G M nxn is a time-invariant system ma¬ 
trix that represents the linear interconnections between agents. Note that 
we have A l} = 0 if there is no interconnection between the agents i and j. 

To represent the entire systems dynamics, we define the state x(k) G 
jgAfnxJVn as A [ a ; 1 (k) T , x 2 (k) T ,.. ., Xn ( k) T ] T . Then, the system dynam¬ 
ics of the DNCS is given as 


x(k + 1) = Ax(k), 

where the matrix A G ]^ NnxNn is defined by 
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0, if no connection between the agents i and j, 
Aij, otherwise. 
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Figure 1 : The schematic of the large-scale distributed networked control system. The 
communication delay is represented by k* = k — t, where r is the random communication 
delay term. 

For the discrete-time system in (|2]), it is well known that the system is 
stable if and only if the condition p(A) < 1 is satisfied. We assume that 
the system in ([2]), which is the case without communication delays is stable 
throughout the paper. Then, we address the problem to analyze the stability 
in the presence of random communication delays. We remind the reader that 
N is very large. 

2.2. DNCS with communication delays 

Often, network communication between agents encounter time delays or 
packet losses while sending and receiving data as described in Fig. [l] We 
denote the symbol r as communication delays and assume that r has a 
discrete value bounded by 0 < r < r d < oo, where Td is a finite-valued 
maximum delay. Then, the dynamics for the agent i with communication 
delays can be expressed as: 

Xj(k + 1 ) = ^ AjjXjjk*), i = 1 , 2 , ( 3 ) 

jeM, 

where k* = k — r. Note that we have no communication delays when i — j 
since there is no communication in this case. 

The communication delay, modeled as a stochastic process, is represented 
by the term k*. To analyze the stability of the DNCS, we define an augmented 
state X(k) as X(k) = [x(k ) T , x(k — 1) T , • • • , x(k — T d ) T ] T e ygNnqxNnq^ w i iere 
q = r d + 1. Then, the dynamics for the entire system is given by 

X(k + 1) = W(k)X(k), (4) 
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where W(k) = 


Mk) 

A-2 (k) 

••• (k) 

A q (k) 


I 

0 

0 

0 


0 

I 

0 

0 

G ]^NnqxNnq 

0 

0 

I 

0 



the matrix / denotes an identity matrix with proper dimensions, and the 
time-varying matrices Aj(k ) G M. NnxNn ; j — 1^2,... ,q, model the random¬ 
ness in the communication delays between neighboring agents. 


3. Switched System Approach 

Without loss of generality, the dynamics of the large-scale DNCS with 
communication delays in Q can be transformed into a switched system 
framework as : 

x(k + 1) = W^ k )x(k), a(k) G {1,2, • • • ,m}, (5) 

where the set of matrices {Wa(i)}™/u =1 represents all possible communica¬ 
tions delays between interconnected agents, {cr(fc)} is the switching sequence, 
and m is the total number of switching modes. When the switching sequence 
(cr(fc)} is stochastic, (J5]) is referred to as a stochastic switched linear system 
or a stochastic jump linear system [7j. For the stochastic switched linear 
system, the switching sequence {cr(fc)} is governed by the mode-occupation 
switching probability n (k) = [7Ti(A:), vr 2 (fc),..., 7r m (k)\, where 7r* is a frac¬ 
tion number, representing the modal probability such that Y^=i n i = 1 and 
0 < 7Tj < 1, Vi. Typically, randomness in communication delays or commu¬ 
nication losses has been modeled by the MJLS framework 0, ra, mi, mi, 
|l 3 j . Therefore, we make the following assumption in our analysis. 

• Assumption: Consider the stochastic jump linear system (J5]) with the 
switching probability 7r (k) = [7r 1 (/c), 7r 2 (fc),..., 7T m (k)\. Then, ir(k) is 
updated by the Markovian process given by n(k + 1) = 7 r(k)P, where 
P g M mxm is the Markov transition probability matrix. 

Since the MJLS is a family of the stochastic switched linear system, various 
stability notions can be defined [H] . In this paper, we will consider the mean 
square stability condition, defined below. 

Definition 3.1. (Definition 1.1 in The MJLS is said to be mean square 
stable if for any initial condition xq and arbitrary initial probability distribu¬ 
tion 7r(0) ; lim E [| \x(k, x 0 )| | 2 ] = 0. 
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Note that for the large-scale DNCS, the total number of switching modes 
m depends on the size q and N. Since the communication delays take place 
independently while receiving and sending the data for each agent, m is 
calculated by counting all possible scenarios to distribute every matrices 
Aij G M nxn for i ^ j in the block matrix A G ^ NnxNn given in ([2]), into each 
Aj(k) G M. NnxNn , j — f ; 2, ..., q, given in Q, which results in m = q N(N k. 
For large N, m is quite large, which makes current analysis tools for the 
MJLS computationally intractable. 

Before we further proceed, we introduce the following proposition that 
was developed for the stability analysis of the MJLS. 


Proposition 3.1. (Theorem 1 in }21f) The MJLS with the Markov transi¬ 
tion probability matrix P is mean square stable if and only if 


p ((P T ® I) diag (Wj ® Wj)) <1, (6) 

where I is an identity matrix with a proper dimension, 


diag (Wj ®Wj) = 


(Wi ® Wi) 
0 


0 

0 


0 

(W 2 ® W 2 ) 

0 

0 


0 

0 


0 

0 


(Bm-l ® W m . l) 0 

0 (W m ® W m ) 


and m is the total number of the switching modes. 


For the given set of matrices {W a ^ fc)}£L.) =1 and the transition probability 
matrix P, one can always compute the spectral radius given in (|6j) , and hence 
guarantee the system stability. 

Unfortunately, this condition is not applicable to large-scale DNCSs since 
N is very high and results in extremely large m. For example, even if q — 2 
and N = 100, we have m = 2 100x99 . It is not possible to compute the 
spectral radius for such problems. To circumvent this scalability issue, we 
present next a new analysis approach for such large-scale DNCSs. 


4. Stability with Reduced Mode Dynamics 

In this section, we define a new augmented state to reduce the mode 
numbers as follows: 

Xifk) = [xi(k) T , Xifk - 1) T , • • • , Xi(k - r d ) T ] T G M nin9 , 
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where x^k) = [xi(k ) T , Xj(fc) T ] T G M nin , h* = and G M n with 

j G A/) denotes all states that are neighbor to xfk) G M n . 

Then, we can construct a switched linear system framework similarly to 
(|5]) as follows: 


Xi(k 

+ 1) = 

Wf i{k) Xi(k), 

(?i{k) G 

{1,2,.. 



'ii(fc) 

A 2 {k) ••• 

A q -i(k) 

A q (k) 



I 

0 

0 

0 


where Wf {k) ± 

0 

/ 

0 

0 



0 

0 

/ 

0 



with the time-varying matrix Aj{k ) G M n ; nxn i n ; j = 1, 2 ,... ,q. In this case, 
the total number of the switching modes for (J7|) is given by = q n A n i~ 1 ) m 
By implementing the reduce mode model given in (J7]), we will provide a 
computationally efficient tool for the stability analysis of the original DNCS 
in the following theorem. 


Theorem 4.1. Consider the large-scale DNCS (J5]) with Markovian commu¬ 
nication delays accompanied by the transition probability matrix P. The nec¬ 
essary and sufficient condition for the mean square stability of this system is 
then given by 

p((P ?;T <g> /)diag(Wj 0 IT])) <1, Vi = 1, 2 ,..., N, (8) 


where P l G is the transition probability matrix for the reduced mode 

MJLS given in (J7]) ; I is an identity matrix with a proper dimension, N is 
the total number of the agents in the system, rn t = q n A ni ~ 1 ) i s the total mode 
numbers for the reduce mode MJLS, and 


'(Wi 0 Wf) 
0 


diagfiVj 0 Wf) A 


0 

0 


0 

(Wl 0 Wf) 

0 

0 


0 

0 


0 

0 


«-i0<-i) 0 

o (wi 0 wf 


Proof. Let the matrix Q l {k) be of the form Q'fk ) = K[xi(k)xi(k) T ], Then, 

mi 

Q l {k ) is alternatively obtained by the following equation: Q l (k ) = E Qi(k), 

s =1 
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where Q l s {k ) = E Xi(k)xi(k) T \<7i(k) = s it l s (k), and vr*(/c) = Pr(<Tj(/c) = s). 
Then, Q\(k) satisfies 


mt 


Ql(k) = ^ ~2E[xi(k)xi(k) T | Ui(k) = s,a.i(k - 1 ) = r] 

Pr(<Ji(k - 1) = r | Gi(k) = s)Tt l s {k) 

mi 

= ^E [xi(k)xi(k) T | Gi{k) = s, ai(k - 1) = r] 

Pr(cq(fc) = s | Gi(k — 1) = r) 7r*(h — 1) 


r= 1 


~Pr 


mi 


\~2pIs E [xi{k)xi(k) T I (Jj(h) = s,Gi(k - 1) = r}P r (k - 1) 


r= 1 
mi 


= ^2p l rs E[W^. (jk _i)Xi(fc - l)xi(fc - l) T Wj. (fc _ 1) I - 1) = r]<(fc - 1) 

r =1 
mi 

= p l rs W l r E [£i(k — 1 )xi{k - 1) T | Gi{k - 1) = r\n l r (k — 1) wf 

< -----' 

=Q*(fc-l) 

mi 

= y2Prs^Qr(k-^)Wf- 

r= 1 

In the second equality of above equation, denotes the mode transition 
probability from r to s in the Markov transition probability matrix P l . 
Taking the vectorization in above equation results in 

mi 


vec (Ql(k)) = vec ^p l rs W l r Q l r {k - 1)W* 

\r =1 

mi 

= 5>t s vec (w'QKk - l)Ihf 

r=l 
mi 

= ® fT>ec(Q;(£; - 1)). 


r=l 


In the last equality, we used the property that vec(ABC) = (C T (8)kL)vec(5). 



We define a new variable y l ^(k) = vec {Q l ^(k)j, which leads to 


yi(k) = '£p < r,(K®w;)yUk- i). 

r —1 

By stacking y l ^{k) from 1 up to m;, with a new definition for the aug¬ 
mented state y l {k ) = [y\(k) T y l 2 (k) T ... y l mt (k ) T ] T , we have the following 
recursion equation: 


y\k) = 


PlAwl^wl) p^iW^Wi) 
p\ 2 {Wi®W{) pAM^M) 


PinAWL^W^) 

PAmAW'^WL). 



~ y\(k-l)~ 


yi(k- 1) 


Vim (k - 1)_ 


= (P iT ®/)diag(W*0W; 


=y i (fe—1) 


From the above equation, it is clear that py(P l <S> J)diag(hFj <8> Wj ) 
1 implies lim = 0, and hence this leads to lim Q l {k) = 0 - 

k—yoG k —^oo 

lim tr ( Q\k )) = 0 4=^ lim E [||£j(fc)|| 2 l = 0, which is the sufficient 

fc—>oo v ' fc—>■ oo L J 


< 


mean 


i.T 


square stability condition for Xi(k). On the other hand, if we have py[P l 

/)diag(Wj ® Wj) j > 1, then y l (k ) will diverge, resulting in necessity for the 

mean square stability of Xi(k). Hence, the spectral radius being less than 
one is the necessary and sufficient mean square stability condition for the 
state Xi(k). Further, we have lim E i||:rj(fc)|| 2 l = 0, Vi = 1,2,..., A/ - 

lim E [||a;(A;)|| 2 ] = 0, where x(k) is the state for the DNCS defined in (|5]). 

fc—»oo L J —' 

This concludes the proof. 


Remark 4.1. Theorem fA. provides an efficient way to analyze the stability 
for the large-scale DNCSs. The key idea stems from the hypothesis that the 
stability of each subsystem by partitioning the original system will provide 
the stability of the entire system. Without any relaxation or conservatism, 
theorem \f.l\ proved the necessary and sufficient condition for stability, which 
is equivalent to © for the mean square stability of the entire system. Com¬ 
pared to the total number of modes of full state model ([5]), which is q N W~ 1 ) > 
the reduced mode model ([7]) has total YliLi q n d n i~T modes. Consequently, the 
growth of mode numbers in full state model is exponential with respect to N 2 , 
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whereas that in reduced mode model is linear with regard to N. Therefore, 
theorem EH is computationally more efficient. 


5. Stability Region and Stability Bound for Uncertain Markov 

Transition Probability Matrix 

The Markov transition probability matrix can be obtained from data of 
communication delays. However, the statistics itself contains uncertainties 
due to the uncertainty in the data. Thus, one cannot estimate the exact 
transition probability matrix in practice. In this subsection, we assume that 
the Markov transition probability matrix has uncertainty, i.e. P l = P*+AP®, 
where P l is the nominal value and A P l is the uncertainty in the Markov 
transition probability matrix for i th subsystem. Due to the variation in A P l , 
the system stability may change and hence we want to estimate the bound 
for A P l to guarantee the system stability. Here we assume that A P l has the 
following structure: 



Since we have a constraint such that the row sum has to be zero for A P l 
in above equation, we aim to find the feasible maximum bound for each row, 
e®, satisfying the inequality |Ap* s | < e®, Vr = 1,2,..., m*, to guarantee the 
system stability. Then, e\ , Vr = 1,2,..., m* can be obtained by the following 
two steps. 

Step 1: Solve via Linear Programming (LP) 


maximize 1 T z (for upper bound) 


( 10 ) 


or minimize 1 T z (for lower bound) J 

subject to A|z| < b s , Vs = 1, 2,..., m* 


lh s < z s < ub s ,Vs = 1,2 ,mi 
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( 11 ) 

( 12 ) 





where 


%a - l A Pls > A P 


2s 5 


, A p. 


[ Z 1 ) ^2 ) 


,* 1 T 

m^sJ ? 


„T iT 


A A 


^1} ^2? ? ^r; 


, with a,- A || W* 0 W■ Hoc, j = 1 , 2 ,..., m h 


(-» = i-E 




r=l 




lbs = [-Pn, -P l 2s, 
ub s A [l - p[ s , l - ... l - 


The inequality constraint in the LP problem guarantees the mean 
square stability according to the Lemma 5.1 and Theorem 5.1| The term 
lb s and ub s in (12) are the lower and upper bounds for z s , according to 


0 < pis = (pL + a pU) < 1- 


Step 2: Obtain Feasible Solution with Hyperplane Constraint 

We can compute the feasible maximum bound for A p l rs as follows. 


£ r = nun 


(min(|4ib|),min( |4 iUb |)), r = 1 , 2 ,...,m,. 


(13) 


where e l r lh A 


4,ub - [( A Prl)ubi ( A PU) ub. 


- [( A 4i)fb> ( A PU)mi • • • 1 (ApiJmVi 

■■■1 ( A 4m) ub] T i and ( A Pis)lb, ( A Ks)lb denote optimal lower and upper 
bounds for A p l rs , obtained from the LP, respectively. 

Since upper or lower bounds are solved by maximizing or minimizing the 
objective function, (A p l rs )* has different values for upper and lower bounds. 
Fig|2] shows the geometry of stability region analysis for uncertain transition 
probability matrix. The region Si stands for the bounds that come from 
— p l rs < A p l rs < 1 — p l rs . S *2 can be obtained from inequality constraint 
0. The region S 3 denotes the solution from the LP and S is the feasible 
maximum upper and lower bounds with a stability guarantee. Note that 
A P l satisfies A Pis = 0, Vr = 1, 2,..., m; and hence feasible solutions 

should lie on the hyperplane, satisfying A p 2 rl + A p l r2 + ... + A p l rrn . = 0, 
Vr = 1, 2,..., raj. Therefore, we can compute the feasible maximum bound 


from (13) for each row r. 


Now we prove that inequality constraint 0 guarantees the system sta¬ 
bility. 
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Figure 2 : The geometry of the Stability Region Analysis for the uncertain Markov transi¬ 
tion probability matrix when nii = 2 . Each region is described in the figure. 

Lemma 5.1. Consider a block matrix X defined by 


Xu 

Ah 2 • 

X im 

X 21 

a 22 • 

^2m 

^ml 

A m 2 

Y 

yv mm 


where matrix e M nxn . Then, we have p{X) < 1, if E IA«IL < i. v< = 

J=1 

1 ,2 ,... ,m. 

PROOF. For the block matrix X given above, the following inequality condi¬ 
tion || X ||oo< max* || Xij Hoc holds. Also, it is well known that p(X) < 

|| X ||p for any choice of p. 

Therefore, we conclude that Xq=i II -^-ij lloo< 1, Vi = 1,2=k 
P(X) <|| X ||oc< 1. 
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Theorem 5.1. Consider the MJLS for the large-scale DNCS with com¬ 
munication delays. Then, ([5]) is mean square stable if 


mi 

^2 a r\&Prs\ < Ps, 
r= 1 


Vs 

Vi 


1.2.. .., nii, 

1.2.. .. ,N 


mi 

where a r =|| If]? 0 IT] ||oo and (3 S = 1 — || IT] 0 IT] Hoc, is satisfied. 

r =1 

PROOF. If the Markov transition probability matrix for the system in ([7]) 
has the uncertainty denoted by P l = P l + A P\ then the term p^(P* T <g) I) 

diag(IT] 0 Wj) j in ((8|) can be expressed as 

p((P iT 0/)diag(IT]0 IT])) 

=p(((P* + AP*) T 0 7)diag(IT] 0 IT])) 

=p(^((P iT 0 /) + (AP* T 0 /))diag(IT] 0 IT])) 

=p(jP* T 0 /)diag(IT] 0 IT]) + (AP' ;T 0 /)diag(IT] 0 IT])) 

< || (P' T 0 /)diag(IT] 0 IT]) + (AP jT 0 J)diag(IT] 0 IT]) IU 

< || (P' T 0 /)diag(IT] 0 Wj) ||oo + || (AP* T 0 J)diag(IT] 0 IT]) (U, (14) 

In the first inequality, we used the fact that p(-) <|| • ||oo and the sub- 
multiplicative property was applied in the last inequality. The block matrix 
structure for each term of the last inequality is alternatively expressed as 
follows: 
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(P ?;T <g> /)diag(Wj <g> Wj) 


' Pi1 1 

P 21 1 • 

■ Pimi 1 " 


~Wf < 8 ) Wf 

0 

0 

PA 

P22 1 ■ 

Pmf 


0 

if ® wi 

0 

firm 1 

P2mj ■ 

Prn.rn, P 


0 

0 

Wf ® 


={P i <g>/) 


=diag(W;<g>W 7 ) 


7iPn 72^21 


IrnfmA 


71P12 72P22 ••• 7m i pA 2 


71 Plrm l-i.p2.rm ■ ■ ■ Irnfrmrm 

OO 



where 7 j = (Wj 0 W]), j = 

1,2,.. 

. , 777 *, and similarly, 




7iAp l n 

72 a jy 21 • • 

• 7m, Apj^! 

T * 


7iAp { 2 

72 Ap 22 • • 

• ImAPma 

(A P l 0 /)diag(lTj 0 Lbj) 

OO 




7iApi m . 

72Ap 2m . • • 

■ ImAplmrrn 


By applying the result in Lemma 5.1 into (|14() , it is guaranteed that 

T 


p( (P l 0 J)diag(hhj 0 1-Lj) J < 1, if the following condition 

mi mi 

J^a r |A/4 + ^a r p rs <1, Vs = 1,2, ...,m;, 


r =1 


r=l 


where a r A | If™ <g) 11 ^, is satisfied. 

Therefore, (J5]) is mean square stable by Theorem 4.1 if it is guaranteed 
that 


m; 


^a r |Aj4| < p a , 


r=1 


VS = 1,2, . . . ,77lj, 

Vi = 1,2,..., iV, 


where = Yfli a rPr 
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6. Examples 

6.1. Stability Analysis for N Inverted Pendulum System 

Consider N inverted pendulum system, which are physically intercon¬ 
nected by linear springs |22j . The discrete-time subsystem dynamics with 
communication delays is modeled by 

Xi{k + 1 ) = AiXi(k) + BiUi(k) + ^ (/,•'). 

jeN* 

with subsystem matrices: 


0 O' 

t oj ’ 

where k denotes the discrete-time index and Xi = (xi x , Xi 2 ) T G M 2 . The 
communication delay is described by the term k* = k — r with the discrete 
value r. The meaning of each parameter and its value are given in Table [lj 
For this system, we consider a state feedback law given by Ui(k) = 
KiXi(k), where K { = [a*/! — rr ^-(8 + -^), —3 ml 2 ] for the control input 

Ui(k*). 


A = 


1 At 
1 


Bi = 


0 
At 
ml 2 


Ha = 


Table 1 : Nomenclature for N Inverted Pendulum Dynamics 


Definition 

Symbol 

Value 

Number of Springs connected to 

di 

1, i = 1,N 

i th Pendulum 

2, otherwise 

Interaction term with neighbours 

hij 

0.04, Vj G Mi 

Gravity 

9 

9.8 

Spring Constant 

K 

5 

Pendulum Mass 

771 

0.5 

Pendulum Length 

i 

1 

Sampling time for discrete-time dynamics 

At 

0.1 


We can rewrite the closed-loop dynamics for this N inverted pendulum 
system as follows: 


Xi(k + 1) = AijXj(k*), where Aj — 
jeATi 


Ai + BiKi, if j = i, 
Hij , otherwise. 
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If there is no communication delay (i.e., k* = k), the dynamics for the 
entire DNCS is given by Q, where the matrix A of which structure is also 
given in Q satisfies p(A) = 0.9525 < 1. Therefore, we can assure that the 
N inverted pendulum system with no communication delays is stable. 

Next, we test the stability for this system with random communication 
delays. We assume that the communication delay r is bounded by 0 < r < 
Td = 1 , i.e., k* = {k, k — 1}, Vi = 1, 2, ..., N. Also, we assume that every 
communication delays are governed by the Markov process with an initial 
probability distribution 7r(0) and the Markov transition probability matrix 
P given by 


7T(0) = [1, 0], P 


0.5 0.5 
0.3 0.7 


(15) 


For this system, even with N = 100, the full state model ([5]) has total 
qN(N-i) _ 2ioox99 rn(X ] es Since this inverted pendulum system has only in¬ 
terconnected terms with neighbors when j = i± 1, otherwise we have Aij = 0. 
Based on this fact and by excluding these cases (i.e., where A t j = 0), we can 
further reduce the mode number to q 2( ~ N - 1 ) — 2 2x ", which is still large. It is 
computationally intractable to deal with 2 2x " numbers of matrices to ana¬ 
lyze system stability. However, in contrast, the reduce mode model ([7]) has 
total YliLi q ni< ' ni ~ 1 ' > = 98 x (2 3x2 ) + 2 x (2 2xl ) = 6280 modes. Furthermore, 
the proposed method fully maximizes its own advantage to reduce the mode 
numbers by considering the symmetric property between agents, which can¬ 
not be implemented on the full state model. Since subsystems are symmetric 
for Vi = 2, 3 ,... , N — 1 and for Vi = 1, N, we only need to check the sta¬ 
bility condition for these two cases. Taking into account the interconnection 
link (i.e., the case where A t] ^ 0), the symmetric structure results in total 
2 2x (3— i) + 2 2x ( 2 ~d — 20 modes, which is drastically reduced when compared 
to 2 2x " numbers of modes. 

The spectral radius for i = 2,3,..., 99 is computed by p[(P l 0 /) 


diag(IF* 0 IF*)) = 0.8864 < 1, where P 1 = (P 0 P 0 P 0 P). For i = 1 and 


N, we have p(jP iT 0 /)diag(IF] 0 IF])J = 0.8682 < 1, where P i = (P0 P). 
Consequently, the N inverted pendulum system is stable in the mean square 
sense according to Theorem |4.1| The state trajectory plot also supports this 


result, as shown in Fig. [3] For this simulation, initial condition was assumed 
to be uniformly distributed in [—1,1], and was generated by manipulating 


16 





>< 

0 

-I—* 

CO 

-I—> 

C/D 



0 10 20 30 40 50 60 70 80 

time(k) 


xT 

0 

-i —' 

CO 

-I—> 

c n 



30 40 50 

time(k) 


Figure 3: State trajectories of each agent for the N inverted pendulum system with the 
Markovian communication delays. Initial conditions are randomly generated for all states 
: Vi, Vi = 1,2,..., 100. 


( Rj 

the MATLAB W command rand(. . .) that generates uniformly distributed 
pseudo random numbers between 0 and 1. 

6.2. Stability Bound for Uncertain Markov transition probability matrix 
In order to solve the LP to estimate the bound for uncertain Markov 
transition probability matrix, we used MATLAB W with CVX[22j, a Matlab- 
based software for convex optimization. 

6.2.1. Scalar system 

Although the proposed method to estimate maximum bound for uncer¬ 
tain Markov transition probability matrix is developed for the large-scale 
DNCS, it is also applicable to general MJLS. We adopted a following ex¬ 
ample, introduced in [IT/] to compare the performance of maximum bound 
estimation. 
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Consider the following M JLS that has two modes with scalar discrete-time 
dynamics. 

x(k + 1) = a a ( k )x(k), cr(fc) G {1,2}, 

a, = 1/2, a 2 = 5/4. 

The Markov transition probability matrix has the form of P = P + A P, 
where 

2 

,J^Ap ij = 0, Vi = 1,2. 

3 = 1 

After applying the two steps proposed in this paper, we obtained the 
maximum bound eq = 0.4, e 2 = 0.02 whereas pT9] gives the value as eq = £ 2 = 
0.021, which is more conservative. For stability check, among all possible 
scenarios with |Ap rs | < e r , \/r,s = 1,2, we have maxp^(P T <g) J)diag(ay $3 

a,j) j = 1, which is a marginal value for stability. Hence, the system is stable 
with obtained maximum bound that is more relaxed than [T9] . 


P = 


0.4 

0.5 


0.6 

0.5 


, A P = 


Apn Ap 12 
Ap 21 Ap 22 


6.2.2. The N Inverted Pendulum System 

Recalling the N inverted pendulum system, we assume that the Markov 
transition probability matrix P* has uncertainty A P l that satisfies AP ? = 
P l — P l . The nominal matrix P l is given by P 1 = (P 0 P) for i — 1, N and 
P l = (P <g) P ® P <g) P) for i — 2, 3,..., N — 1, where P has a same structure 
with the transition probability matrix given in (15). 

The feasible solution with the LP provides the maximum bound e l = 
[4, 4, ...e\ 6 ] = 10" 3 x [4.9, 0.9, 0.9, 0.8, 0.9, 0.8, 0.8, 6.9, 0.9, 0.8, 0.8, 
6.9, 0.8, 6.9, 6.9, 13.5], V? = 2, 3,..., N — 1. For i = 1 and A, we obtained 
e l = [ 4 , 4, 4, 4] = [0-01, 0.01, 0.01, 0.03]. Therefore, we can assure 
that N inverted pendulum system is mean square stable if the uncertainty 
in the Markov transition probability matrix is within the bound such that 
|Ap} s | < er®, Vr and Vi = 1, 2, ..., N. 


7. Conclusions 

This paper studied the mean square stability of the large-scale DNCSs. 
Since the number of modes in such systems is extremely large, current stabil¬ 
ity analysis tools are intractable. To avoid this scalability problem, we pro¬ 
vided a new analysis framework, which incorporates a reduced mode model 
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that scales linearly with respect to the number of subsystems. Additionally, 
we presented a new method to estimate bounds for uncertain Markov transi¬ 
tion probability matrix for which system stability is guaranteed. We showed 
that this method is less conservative than those proposed in the literature. 
The validity of the proposed methods were verified using an example based 
on interconnected inverted pendulums. 
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